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Introduction: Computerized simulations of cardiac activity have significantly contributed 
to our understanding of cardiac electrophysiology, but techniques of simulations based 
on patient-acquired data remain in their infancy. We sought to integrate data acquired 
from human electrophysiological studies into patient-specific models, and validated this 
approach by testing whether electrophysiological responses to sequential premature 
stimuli could be predicted in a quantitatively accurate manner. 

Methods: Eleven patients with structurally normal hearts underwent electrophysiological 
studies. Semi-automated analysis was used to reconstruct activation and repolarization 
dynamics for each electrode. This S2 extrastimuli data was used to inform individualized 
models of cardiac conduction, including a novel derivation of conduction velocity 
restitution. Activation dynamics of multiple premature extrastimuli were then predicted 
from this model and compared against measured patient data as well as data derived 
from the ten-Tusscher cell-ionic model. 

Results: Activation dynamics following a premature S3 were significantly different from 
those after an S2. Patient specific models demonstrated accurate prediction of the S3 
activation wave, (Pearson's R 2 = 0.90, median error 4%). Examination of the modeled 
conduction dynamics allowed inferences into the spatial dispersion of activation delay. 
Further validation was performed against data from the ten-Tusscher cell-ionic model, with 
our model accurately recapitulating predictions of repolarization times (ff 2 = 0.99). 

Conclusions: Simulations based on clinically acquired data can be used to successfully 
predict complex activation patterns following sequential extrastimuli. Such modeling 
techniques may be useful as a method of incorporation of clinical data into predictive 
models. 

Keywords: conduction velocity restitution, computational modeling, action potential duration, patient specific 
modeling, cardiac arrhythmia 



INTRODUCTION 

Computerized simulations of cardiac activity have significantly 
contributed to our understanding of cardiac electrophysiology 
(Carusi et al, 2012). The recent leaps in raw processing power and 
access to supercomputing by research organizations has allowed 
proof-in-principle of many theories of cardiac conduction, from 
the cellular action potential to the generation of arrhythmia in 
heterogenous systems (e.g., in inherited arrhythmia syndromes 
such as Brugada or LQT syndromes). Though simulation has 
been very successful at bridging a knowledge gap between basic 
research findings and understanding of arrhythmia, efforts have 
principally concentrated on the "forward" solution, i.e., pro- 
duction of computer models that inform us of the behavior 
expected from our mathematical knowledge of cardiac physiol- 
ogy (Trayanova, 2011; Krummen et al, 2012). Few studies have 

Abbreviations: ARI, Activation Recovery Interval; APD, Action potential dura- 
tion; CV, Conduction Velocity; AT, Activation Time; ECG, Electrocardiogram; SI, 
S2, S3, Stimulus, Extrastimulus, Second Extrastimulus; ERP, Effective Refractory 
Period; RT, Repolarization Time; MRI, Magnetic Resonance Imaging; VT, 
Ventricular Tachycardia; VF, Ventricular Fibrillation. 



attempted to broach a reverse solution, i.e., the construction of 
an accurate simulation based on patient-acquired data. There 
have been attempts to use patient anatomical data, such as those 
acquired from cardiac MRI, to fit a pre-existing cellular model, 
and some early work has attempted to fit clinically acquired data 
to a simulation (Relan et al., 2011). There remains some dis- 
tance between computational modeling and clinical appreciation 
of arrhythmia generation, and computational approaches have 
not as yet found a clinical role. 

Recently work has illustrated an approach by which elec- 
trophysiological data may be incorporated into simulations. 
Gilmour's group has concentrated on examination of one- 
dimensional conduction models, and has shown how sequential 
close-coupled activations may result in functional block "at-a- 
distance" from the stimulation site i.e., how an activating wave- 
front may impinge upon a prior wave of activation (Gilmour 
et al., 2007; Otani, 2007). Such functional block and wavebreak 
appear to be necessary events in triggering and sustaining the 
development of chaotic human arrhythmias, particularly VF. 
These frameworks incorporate both action potential duration 
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(APD) restitution and conduction velocity restitution, and are 
founded on the cable theory principles originating in the work 
of Hodgkin and Huxley (Noble, 1962). In brief, this concept 
embodies cardiac conduction as a syncytium; strips of cardiac 
muscle are considered to act as uniform cables, allowing con- 
duction in all directions. Unlike the squid giant axon, the cable 
properties of cardiac cells span many different cells, but they 
may still be represented as a two syncytial domains: intracellu- 
lar and extracellular, and conduction properties depend on the 
properties of the cell membranes, ionic concentrations and cel- 
lular geometric arrangements. This concept can be developed to 
higher dimensions e.g., Spach et al., 1979, and/or simplified into 
descriptions of pure propagation of activation and repolarization 
e.g., Nolasco and Dahlen, 1968. These latter coupled-map models 
are easy to deal with both analytically and numerically and avoid 
the large complexities of high-dimensional, non-linear models of 
cardiac conduction. Further iterations of these models have clari- 
fied mathematical descriptions of conduction block at-a-distance 
in one-dimensional fiber (Fox, 2002; Fox et al., 2003). 

Such a model underwent a qualitative validation in a canine 
model by Gelzer et al. (2008), where it was used to predict the 
initiation of VF by combinations of up to 5 extrastimuli. This 
work was notable in that the computer model had incorporated 
biologically-acquired data (action potential restitution measured 
from the interventricular septum). Thus, it represented a proof- 
of-concept of animal-specific electrophysiological models in the 
prediction of VF induction by extrastimuli, and has begun to 
bridge the gap between theoretical electrophysiological modeling 
and clinical application. But whether such models can be applied 
to the human heart, particularly in the investigation of cardiac 
risk, remains unproven. 

A major limitation to closing the gap between ionic car- 
diac models and clinical implementation is the iterative pro- 
cess required to be undertaken to fit the model to the data 
(Clayton et al, 2011). The great number of variables in car- 
diac ionic models precludes accurately fitting of clinical data 
with confidence (Zaniboni et al., 2010). Previous attempts have 
concentrated on anatomical approaches to fit data to the indi- 
vidual, relying on "generic" models of the electrophysiology 
itself (Trayanova, 2011). Recognizing this, we have endeavored 
to describe observed cardiac electrophysiology at the tissue level 
using simple parameters. 

This paper describes the approach we have taken to derive 
computer models from clinical data in our attempt to close 
the gap from theoretical modeling. The mathematical founda- 
tion of our work remains that of the one-dimensional model, 
which is itself a manifestation of cable theory. Thus, this work 
incorporating observed clinical parameters into a predictive 
model, and the validation of this method, in itself acts as 
a validation-step to the relevance of cable theory to human 
arrhythmogenisis. 

We hypothesized that individual human restitution dynamics 
could be described by a one-dimensional model incorporating 
both conduction velocity restitution and action potential resti- 
tution. We tested this by incorporating patient-acquired data 
into individualized one-dimensional models (Hand and Griffith, 
2010). These models were used to predict activation time (AT) 



dynamics following a second extrastimulus, and predictions were 
compared to experimental results. 

METHODS 

PATIENT POPULATION AND ELECTROPHYSIOLOGICAL STUDIES 

Studies were performed in patients undergoing cardiac electro- 
physiological studies with a view to diagnosis and ablation of 
supraventricular tachycardia. All patients gave prior informed 
consent. Studies were performed under minimal conscious seda- 
tion in the post-absorptive state. Antidysrhythmic drugs were 
discontinued for 5 days prior to the study. The study protocol 
had local ethics committee approval and conformed to the stan- 
dards set by the Declaration of Helsinki. Decapolar catheters were 
placed in an epicardial coronary vein via the coronary sinus, at the 
RV apex and retrogradely within the LV cavity adjacent to the epi- 
cardial catheter (St Jude Pathfinder). A reference anodal electrode 
was placed in the inferior vena cava. Further catheters were placed 
according to clinical requirements (Figure 1). Electrograms were 
digitized and recorded at 1000 Hz (Bard Clearsign, MN, USA). 
If no arrhythmia was initiated by ventricular pacing during the 
initial clinical study, restitution studies were performed prior to 
further clinical testing. 

STUDY PROTOCOL 

Following 3 min of steady-state pacing at the RV apex at a cycle 
length of 600 ms, a drivetrain of 10 beats was followed by an 
extrastimulus (S2). The S1S2 interval was decremented by 20 to 
300 ms, thereafter in 5 ms intervals until ERP was reached. The 
ERP was found by increasing CIs by 8 ms and then decrement- 
ing by further 2 ms intervals. In patients undergoing second 
extrastimuli (S3) studies, the S2 extrastimulus was set at 10ms 




FIGURE 1 I Unipolar EGM analysis. The three final of a restitution train 
are shown, typical unipolar electrograms (solid line) and their differential 
(dashed line). Local activation was defined as the most rapid downstroke 
(dVdtmin, circles, A1, A2) of the activation inscription, repolarization was 
defined using the classical (Wyatt) method as dVdtmax of the T-wave 
(asterisks, R-\ , R2). Activation time was defined as the time from pacing 
stimulus (not shown) to local activation, ARI as time from local activation to 
repolarization. It was not possible to reliably measure the diastolic interval 
for close coupled beats, therefore Dl was taken as A^2 — ARI SS . 
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above the observed refractory period, and an S3 was delivered at 
CIs progressively decrementing to a minimum of 180 ms. 

WAVEFORM ANALYSIS 

All electrogram analysis was performed using custom semi- 
automated software running in Matlab r2009b (The Mathworks 
Inc, MA, USA). Local AT and repolarization time (RT) were 
identified using an algorithm and manually checked. Filter set- 
tings of 0.1-150 Hz band-pass (for AT) and 40 Hz low-pass (RT) 
were used. Any electrogram in which the T-wave morphology was 
either inconsistent throughout the drivetrains or inadequately 
defined was excluded from analysis. 

The electrograms of the ultimate 2 beats of each drive train and 
subsequent S2 were exported and analyzed. It is technically chal- 
lenging to measure action potentials in humans at more than one 
site simultaneously. However, the activation-repolarization index 
(ARI), i.e., the time from local activation to local repolarization as 
calculated from a unipolar electrogram signal, is a well-validated 
surrogate for the APD. Techniques for derivation of ARI analy- 
sis from unipolar electrograms have been well-described (Franz 
et al, 1987; Yue et al., 2004; Hanson et al., 2009). Briefly, AT was 
defined at the steepest negative slope of the activation waveform 
(Figure 1), and local repolarization time defined as the steepest 
positive slope of the T-wave (Figure 1), as supported by in-depth 
analysis (Yue et al, 2004; Yue, 2006; Potse et al, 2007, 2009). The 
diastolic interval (DI) preceding S2 was calculated using measured 
values taken from steady state data. Similarly, the DI prior to a S3 
beat was estimated using a restitution model previously obtained 
for each electrode location. AT and repolarization time (RT) were 
measured from the stimulus time (defined as t = 0). The ARI of 
each beat is defined as the time difference from AT to RT, such 
that the repolarization time RT = AT + ARI. 

ARI was plotted against DI to give standard ARI restitution 
curves. AT was plotted against coupling interval to illustrate 
conduction delay. These graphs were used for model creation. 
All electrogram analysis was performed using custom semi- 
automated software running in Matlab r2009b (The Mathworks 
Inc, MA, USA). AT and RT were identified using an algorithm 
and manually checked. Filter settings of 0.1-150 Hz band-pass 
(for AT) and 40 Hz low-pass (RT) were used. Any electrogram in 
which the T-wave morphology was inconsistent throughout the 
drivetrains or inadequately defined was excluded from analysis. 

ONE-DIMENSIONAL CONDUCTION MODEL AND CV RESTITUTION 
MODELING 

We derived a novel method whereby conduction velocity resti- 
tution can be quantified from biologically acquired data using a 
linear cell fiber model. Arbitrary length of cells of 0.1mm was 
used for model derivation. Each cell within a fiber was assumed 
to have homogenous CV and APD restitution properties. Both of 
these factors were described by exponential functions. 

If the length of a tissue unit (e.g., a myocyte) is assumed as 
Ax, and the conduction velocity within tissue as CV, the AT and 
repolarization time of a cell x, can be calculated as: 

Ax 

Activation time: AT n (x,) = AT n (x,_ 1 ) + ( 1 ) 

CV n (xi) 

Repolarization time: _RT„(x,) = AT„(x,) + ARI„(x,) (2) 



Where: AT and RT are the AT and repolarization time; ARI is the 
activation recovery interval; n is the number of beat. 

Hence, the AT and recovery time at cell x; can be expressed as: 



AT„(x,) = AT„(x 0 ) + 



Ax 



;=0 



CV n (Xj) 



RT n ( Xi ) = AT n (x 0 ) + V — 7 . + ARI n ( Xi ) (4) 

— ' ( .V..( r.-i 



(3) 



;=0 



CV n ( Xj ) 



Where: xq represents the site of stimulation. 

Thus, the period between activation n and activation n + 1 for 
site xi can be calculated as: 



AT n+ i (x,) - AT„(x,) = AT n+ 1 (x 0 ) - AT„+ 1 (x 0 ) 



; = 0 



Ax 



Ax 



CV n+1 ( Xi ) CV n {Xi) 



(5) 



Here, the period between activation n and activation n + 1 for site 
xq equals to the coupling interval (C/,,+1 ) applied at site Xn 



AT n+ 1 (x 0 ) - AT„(xo) = CI n+ 1 



(6) 



The period between activation n and activation n + 1 for site x; 
can also be calculated as: 

AT n + 1 (x,) - AT n (*;) = ARI„ (x,) + DI n + 1 (Xi) (7) 

Combining equation 5, 6, and 7, the basic equation for one- 
dimensional conduction model can be expressed as: 



E' / Ax 
7^7 T 
• =0 \CV n+1 (Xj 



Ax 



) CV n {xj) 



■ARI n (xi) 



(8) 



The conduction velocity CV(x,) is assumed to be a function of its 
preceding DI: 



CV(Xi) = CV SS - B CV ■ e' 



Ccv-DI 



(9) 



where CVss is the steady-state conduction velocity, a constant 
defined to be 0.72 (Fox, 2002); Bcv and Ccv are two constants 
determining the steepness of CV restitution curve. This derivation 
is summarized in Figure 2. 

Referring to the experimental data, DI n+ \ (x,) is the DI of pre- 
mature beat S2, ARI n is the ARI of the steady-state activation 
following Si (ARI n = ARIss), CI n + \ is the S1S2 coupling interval. 
These three parameters were directly obtained from experimen- 
tal data. CV n (Xi) is the conduction velocity of a steady-state 
activation, i.e., is equal to CVss- 

For each investigated site, the activation conduction path- 
way from the pacing site to the investigated site was assumed 
to be an independent one-dimensional conduction pathway 
with uniform restitution properties. From the experimental 
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FIGURE 2 | Geometry relationship based on CV restitution. A geometric 
relationship based on the CV restitution model is shown, x,- is the position 
along a conduction path, AD/ is the change in diastolic interval, At is the 
time that an activation wave takes to travel the distance Ax. AT n : 
steady-state activation time with long Dl; RT n : steady-state repolarization 
time; /A7" n+1 : activation time of an early beat with shorter pacing interval and 
consequently showing some CV decrease (CV n+1 (x,-), but as the distance 
over which this operates (Ax) is very small this can be considered as a 
straight line. 



data, an assumed length of the conduction pathway (given the 
assumed steady state conduction velocity) can be calculated by 
the product of steady-state conduction velocity and the steady- 
state conduction time, AT$s(xi), from the pacing site to the 
investigated site: 

L(x,) = CV SS ■ AT ss (x,) (10) 
Hence, the number of cells included in the conduction pathway is 

■ _ L(Xj) 
1 — Ax ■ 

Equation 8 can be rearranged into the following form, which 
allows its solution using experimental data: 

DI S2 (Xi) +ARI ss (x l ) - CI S2 + AT ss (xi) = V A * (11) 

All the dynamics in the left hand side of equation 11 are 
those obtained from experimental data. For each investi- 
gated site, the one-dimensional conduction model and CV 
restitution model can be fitted with its dynamics from 
experimental data. An iterative approach was taken, adjust- 
ing Bcv and Ccv in equation 9, and a minimum least- 
squares error method was employed to select optimum 
parameters. 

INCORPORATION OF PATIENT DATA INTO ONE DIMENSIONAL MODELS 

Data acquired from patient studies was applied to the model, and 
a separate (and independent) model was created for each elec- 
trode site. ARI restitution was modeled using a mono-exponential 
equation of the form specified in equation 13. 

A standard value of 0.72 ms-1 was assumed as baseline con- 
duction velocity, changing this absolute value by ±0.2 ms did not 
significantly affect derived CV dynamics. Despite model fitting 
being computationally intense, the use of simplified exponen- 
tial dynamics allowed modern desktop machines to derive CV 



parameters reasonably quickly. For a typical patient, 100-2000 
cells were included, depending on the distance of the con- 
duction pathway. For example, if the steady-state conduction 
time was 100 ms at the investigated site, assuming the steady- 
state conduction velocity to be 0.72 ms-1 and cell length to be 
100 |im, the number of cells included in the conduction model 
was 720. 

Validation of our method was performed in two ways, com- 
parison with ionic cellular model system and comparison with 
experimentally acquired patient data. 

Once both ARI and CV restitution parameters had been 
derived from patient data, the cell fiber model was used to predict 
conduction dynamics of a second extrastimuli (S3) introduced 
following a short S1S2 coupling interval. Modeling of sequential 
beats was only attempted if the r 2 for fitting to S2 data was >0.7 
at that electrode site, and sites where electrograms were of poor 
quality were excluded. ATs following the S3 beat obtained experi- 
mentally were qualitatively compared to those obtained from the 
cell fiber model. 

COMPARING THE EXPONENTIAL MODEL OF REPOLARIZATION 
DYNAMICS AGAINST DATA IONIC CELLULAR SIMULATION 

At short S2S3 coupling intervals, the QRS complex of the unipolar 
electrogram following an S3 can widen to fuse with the repolar- 
ization complex. This precludes using unipolar electrograms to 
validate our simple model for prediction of repolarization times. 
We thus performed preliminary validation of our model for pre- 
diction of S3 repolarization times against data taken from an 
ionic cellular model. The ion channel model is obtained from ten 
Tusscher et al. from their studies of human ventricular tissue (ten 
Tusscher et al, 2009). 

The ten Tusscher human ventricular ionic model includes 16 
types of currents, and 30 ionic gate variables. It is practically 
impossible with current technology to obtain detailed informa- 
tion of each type of ion current or gate variables on an individual 
patient. However, the experimental data provides information on 
key cellular properties such as APD, AT, and repolarization time. 
These three interdependent properties were fitted with the resti- 
tution models, which is then used to simulate interactions of 
activation and repolarization. 

The ten Tusscher human ventricular model can be summa- 
rized in a general form below (the forward Euler method was 
used to integrate the model). Initial conditions for ionic values 
and conductances are given in Table 1 : 

oV _ Ii on + Istim 

dt C m 

Where: V is voltage, 
t is time, 

lion is the sum of all transmembrane ionic currents, 

7 stlm is the stimulus current, 

C m is cell capacitance per unit surface area. 

Action potential duration (APD) is defined to be action poten- 
tial duration at 90% repolarization (APD90). 

Diastolic interval (DI) equals to the coupling interval of 
S2 minus APD from steady-state beat (APD 5S ). 
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Table 1 | Initial conditions for state variables in human ventricular ionic model. 
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Where; V, membrane voltage; [/Va+],-, intracellular sodium ion concentration; lNa + ] ss , concentration of sodium ion, in subspace compartment; [K + ]-,, intracellular 
potassium ion concentration; IK + ] SS , concentration of potassium ion, in subspace compartment; [Ca 2 +];, intracellular calcium ion concentration; [Ca 2 +J ss , concen- 
tration of calcium ion, in subspace compartment; [Ca 2+ ]„ sr , concentration of calcium ion, in network SR compartment; [Ca 2 +]p r , concentration of calcium ion, in 
junctional SR compartment; m, activation gate for fast Na + current l^ a ; hf ast , fast development of inactivation gate for fast Na + current ln a ; h s i ow , slow development 
of inactivation gate for fast Na + current l^ a ; j, recovery from inactivation gate for fast Na + current ln a ; hc a MK,siow, slow development of inactivation gate for CaMK 
phosphorylated fast Na + current lnai icaMK, recovery from inactivation gate for CaMK phosphorylated fast Na + current l Na ; m L , activation gate for late Wa+ current 
Iftai l~>L, inactivation gate for late Na + current l^ a ; hic a MK, inactivation gate for CaMK phosphorylated Na + current l^ a ; a, activation gate for transient outward K + 
current l to ; i fast , fast inactivation gate for transient outward K + current l to ; i s i ow , slow inactivation gate for transient outward K + current l t0 ; a^MK, activation gate 
for CaMK phosphorylated transient outward K + current l t0 ; ic a MK.f a st, fast inactivation gate for CaMK phosphorylated transient outward K + current l to ; ic a MK,siow, 
slow inactivation gate for CaMK phosphorylated transient outward K + current l to ; d, activation gate for Ca 1+ current through the L-type Ca 2+ channel lent, ffast, 
fast voltage dependent inactivation gate for Ca 2 ^ current through the L-type Ca 2+ channel l CaL ; f s i ow , slow voltage dependent inactivation gate for Ca 2+ current 
through the L-type Ca 2+ channel l CaL ; fca fasti f ast development of CcP- + dependent inactivation gate for Ca 2+ current through the L-type Ca 2+ channel / Cat ; fca, s iow 
slow development of Ca 2 ^ dependent inactivation gate for Ca 2+ current through the L-type Ca 2+ channel IcuLiica, recovery from Ca 2+ dependent inactivation gate 
for Ca 2+ current through the L-type Ca 2+ channel l Ca C n, fraction in Ca 2+ dependent inactivation mode for Ca 2+ current through the L-type Ca 2+ channel l CaL ; 
fc a MK,f a st, f as t development of Ca 2+ dependent inactivation gate for CaMK phosphorylated l Ca j_; fc a ,c a MK,f a st* slow development of Ca 2+ dependent inactivation 
gate for CaMK phosphorylated lc a C, *r,f a st, fast activation/deactivation gate for rapid delayed rectifier K + current I&; x r s i ow , slow activation/deactivation gate for 
rapid delayed rectifier K + current l Kr ; x sl , activation gate for slow delayed rectifier K + current l Ks ; x S 2, deactivation gate for slow delayed rectifier K + current l Ks ; 
x K i , inactivation gate for inward rectifier K + current fa ; J re / /yp ( non-phosphorylated Ca 1+ release, via ryanodine receptors, from jsr to myoplasm; J re i,c a MK, CaMK 
phosphorylated Ca 2+ release, via ryanodine receptors, from jsr to myoplasm; CaMK trap , fraction of autonomous CaMK binding sites with trapped calmodulin. 



The form of exponential model tested against data acquired 
from the ionic cellular model was: 

APD = APD ss _ Bme c.Di (13) 

The pacing protocols applied to ion channel model were the 
same as those used clinically (i.e., a sequential extrastimulus 
study). Data acquired from the ion channel model was com- 
pared against S3 repolarization times modeled using our simple 
exponential-based approach. 

STATISTICS 

Data is presented as mean ± standard deviation unless oth- 
erwise stated. Comparative statistics were calculated within R. 
An individual simulation was performed for each measured site, 
thus each patient had multiple simulations performed. Multilevel 
regression enabled appreciation of the statistical reliability of pre- 
dictions without inappropriate prejudice or favor regarding the 
grouped and leveled nature of measurements (i.e., incorporating 
both within-subject and between-subject variances). These cal- 
culations were performed using the lmer package (Pinheiro et al, 



2009), single-level comparisons of normally-distributed contin- 
uous variables were performed using Student's T-test. A p-value 
of less than 0.05 was taken as significant. Correlation coefficients 
involving model comparisons were performed with Matlab. 

RESULTS 

PATIENT DEMOGRAPHICS 

Eleven patients were enrolled in the study, aged 29 ± 10 years, 
7 females. All patients had structurally normal hearts on echocar- 
diography, and underwent EP studies for symptoms consistent 
with paroxysmal supraventricular tachycardia. Final diagnoses 
were atrioventricular re-entrant tachycardia (5 patients), parox- 
ysmal atrial fibrillation (1). No arrhythmia was induced in 5 
patients. Six of these patients underwent S3 studies and were 
included in the S3 validation. 

ACTIVATION TIME, ARI RESTITUTION AND REPOLARIZATION CURVES 

The dynamics of how AT, ARI, repolarization time (RT) and 
DI varied with the coupling interval of the first extrastimulus 
were examined. This is illustrated in Figure 3. As coupling inter- 
val decreases from steady state, ARI restitution is observed, and 
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^r J2 (x M ) AT sl (x,) RT s2 (.x,_,lRT s 4x,) 



FIGURE 3 | Geometry relationship and derivation of modulation of 
dispersion. A premature extrastimulus results in conduction velocity 
restitution. Consequently, ATs are delayed, but this delay allows resumption 
toward steady state conduction velocities over the course of the activation 
path. The modulation of activation dispersion is amplified via ARI restitution 
on the repolarization wavefront, and dispersion of repolarization increases 
further. 



both ARI and repolarization time shorten. But as coupling inter- 
vals shortened further toward ERP, the effects of conduction 
velocity restitution became apparent (as illustrated by increasing 
ATs). This prevents further decrease in local DIs, and hence ARI 
restitution is restricted. The local repolarization time effectively 
represents the cumulative effect of ARI and conduction velocity 
restitution, together these produced a distinctive repolarization 
time curve illustrated in Figure 4. 

MODEL FITTING TO S 2 DATA 

One-dimensional models were fitted to clinically acquired S2 
data in all 11 patients. ARI restitution curves were accurately 
described by a simple exponential model (Coefficient of determi- 
nation, r 2 = 0.966 [95% CI 0.963-0.968], p < 0.0001). AT curves 
were similarly accurately reproduced following our derivations 
of conduction velocity restitution, (r 2 = 0.993 [0.992-0.993], 
p < 0.0001). 

ACTIVATION TIMES FOLLOWING SEQUENTIAL PREMATURE STIMULI 

To show that prediction of ATs following sequential extrastimuli 
was a reasonable standard by which to assess the accuracy of our 
model, we compared the activation dynamics following S2 and 
S3 stimuli. At similar coupling intervals, AT following a prema- 
ture S3 was significantly shorter than after an S2 [AT normalized 
to steady state: 1.1 following S3 (95% confidence interval 1.1- 
1.2) vs. 1.6 following S 2 (1.5-1.7), p < 0.0001]. ERP of S 3 was 
shorter than ERP of S 2 (<186 ± 8 vs. 222 ± 25 ms, p < 0.0001). 
More activation delay existed after an S3 pre-ERP than after an S2 
(39 ± 28 vs. 33 ± 16ms,p = 0.0003). Activation dynamics were 
thus significantly different following S2 and S3 stimuli, with the 
ERP of an S2 considerably shorter than that following the steady 
state beats. 

VALIDATION OF ID MODELS 

We tested whether the one-dimensional model, incorporating 
both ARI and CV restitution, could explain and predict the 
observed variation in activation dynamics following S2 and S3 
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FIGURE 4 | Interactions of activation and ARI leading to repolarization 
time change. The response of diastolic interval (Dl), activation time and 
ARI to shortening coupling interval are shown, with the overall 
repolarization time dynamic shown in the top panel. As paced intervals 
decrease, diastolic intervals decrease initially in a linear fashion (Bottom 
panel). The reduction in Dl is associated with a non-linear reduction in ARI 
(restitution). As coupling intervals become very short, activation time 
increases due to the effect of conduction velocity restitution acting 
between the stimulus site and the recording electrode. Activation time 
increases prevent further reduction of local diastolic interval (Bottom 
panel), which is seen to level off. ARI restitution is blunted at these short 
coupling intervals. Repolarization time thus increases as a result of 
increasing activation times. Representative patient data from RV 
(Star — early activating site, cross- mid activating site, circle- late activating 
site, dotted lines represent smoothed dynamics). 



stimuli. This serves as a quantative validation of the ID model 
incorporating CV and ARI restitution. Validation was performed 
using data from 6 patients who underwent S3 protocol. The 
coefficient of determination (predicted values vs. measured), R 2 , 
was 0.902 [95% CI: 0.890-0.911] (p < 0.0001, Figure5). The 
high predictive accuracy possible with this method serves as a 
quantitative validation of this one-dimensional model of cardiac 
activation in human hearts. 

SIMULATING MODULATION OF DISPERSION OF REPOLARIZATION 

Simulations based on patient data enable examination of elec- 
trophysiology which otherwise remains impenetrable through 
traditional experimental methods. We examined how dispersion 
of repolarization varied with shortening of S 2 coupling interval. 
True dispersion of repolarization, the difference from the earliest 
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FIGURE 5 | Correlation of measured and predicted AT of second 
extrastimulus. The correlation of modeled S3 activation time with 
measured S3 activation is shown. Each patient is represented by individual 
symbols, r 2 = 0.90. 



repolarization to the latest repolarization, has both a spatial and 
temporal component. It is very hard to measure repolarization 
time accurately very close to the stimulus site due to stimulus 
artifact. But the most significant repolarization gradients will be 
close to the stimulus site following a close-coupled S2; the area 
immediately at the stimulus site is exposed to very little delay 
prior to activation, whereas further sites only a very short distance 
away experience significant activation delay. A patient-based sim- 
ulation allows appreciation of these gradients (Figure 6). In the 
example shown the dispersion of repolarization increased from 
41 to 70 ms between mid to late activating sites, yet increased 4- 
fold from 25 to 100 ms between earlier activating sites. Thus, the 
induction of repolarization gradients by premature stimuli occurs 
very close to the stimulus site. 

VALIDATION OF REPOLARIZATION TIME PREDICTIONS AGAINST IONIC 
CELLULAR SIMULATIONS 

For S3 responses, there was close agreement between repolariza- 
tion times determined from our model and repolarization times 
acquired using the ionic model (coefficient of determination, 
r 2 = 0.997 [95% CI 0.993-0.999], p < 0.0001). This supports the 
validity of our model in comparison to full ionic cellular models 
undergoing sequential premature extrastimuli. 

DISCUSSION 

Using a combination of physiological measurement and simula- 
tion, this study demonstrates that the complex dynamic interac- 
tions between activation and repolarization induced by sequential 
extrastimuli can approximated in terms of ARI (APD) and CV 
restitution. This serves as a strong validation of a personal- 
ized computational approach incorporating biologically acquired 
variables and it serves as an extension and validation in humans 
of previous animal work (Gilmour et al., 2007). Our technique 
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provides a platform by which individual patient data may be 
incorporated into computerized simulations and serves as an 
example of a method by which the gap between fundamental 
cardiac simulation and clinical science can be bridged. 

This study is the first to quantitatively validate a one- 
dimensional model of cardiac conduction in man. Such patient- 
based models provide a method by which physiologically 
acquired data may be interrogated in a computerized simulation. 
We have shown that endocardial physiology can be adequately 
modeled by simple mathematical descriptors of CV and ARI 
restitution curves, and confirm the relevance of one-dimensional 
conduction models incorporating both these features to human 
cardiac electrophysiology. Our validation of a straightforward 
one-dimension model can allow further in-depth examination 
of human physiology in a manner that would not be possible 
experimentally, for both ethical and practical reasons. 

The presented modeling is founded on the work of Otani, 
Gilmour and Fox (Fox et al., 2003; Gilmour et al., 2007; Otani, 
2007; Gelzer et al., 2008, 2010). Seminal work by this group, along 
with the single fiber models of Cherry and Fenton (2004, 2007), 
provided initial validation of the use of one-dimensional con- 
duction simulation in biological systems. The ability of such a 
model to predict whether large variations of long-short CIs would 
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induce VF in normal and cardiomyopathic canines was a powerful 
demonstration of the physiological relevance of such simulations 
(Gelzer etal.,2008, 2010). 

We have enhanced these models by fully deriving CV resti- 
tution; allowing better approximation to the actual behavior of 
the activation wavefront and observation of the heterogeneity 
introduced by sequential extrastimuli. This adds to our ear- 
lier description of the modulated dispersion of repolarization in 
humans (Hanson et al., 2009), and addresses questions raised 
from this work (Poelzing and Rosenbaum, 2009). CV restitution 
is often assumed to operate over the entire conduction path in 
a similar manner, obscuring the importance of the spatial dif- 
ferences in CV restitution. Our validation also sought a high 
fidelity correlation of ATs rather than a binary outcome. The large 
datasets acquired from the multiple recording electrodes used in 
each patient (rather than only one in the canine studies) enabled 
an accurate validation of this concept in humans. 

The degree of latency has been described as being important in 
determining the initiation of VT or VF due to programmed elec- 
trical stimulation in earlier clinical studies (Avitall et al., 1992), 
yet the mechanism of this has been unclear. Our explanation of 
CV restitution, and of a mechanism by which cumulative inter- 
actions can set up conduction heterogeneities is consistent with 
these previous descriptions. 

Knowledge of these timing dynamics may help determine sus- 
ceptibility to arrhythmia in channelopathies (Lambiase et al, 
2009; Nam et al., 2010) or a suspected cardiomyopathy. The 
demonstration that patient-specific modeling can predict initia- 
tion sites of clinical arrhythmia underlines the importance of ARI 
and CV interactions in human arrhythmogenisis. 

STUDY LIMITATIONS 

We used clinically acquired data to observe and explain the com- 
plex interactions seen following sequential extrastimuli. These are 
considered central to the initiation of cardiac arrhythmia. We 
did not expect to, nor did we induce arrhythmia in any patient 
undergoing prospective evaluation, and the contribution of these 
mechanisms in clinical arrhythmogenisis was not fully examined. 
Prospective patient data was acquired from contact electrodes; 
simultaneous mapping from the entire endo- and epicardial sur- 
face is unfeasible in conscious patients. However, as a result of 
lack of 3D geometry, conduction path distances were derived from 
timing data and knowledge of electrode spacings. The demonstra- 
tion that these techniques can effectively predict wavebreak and 
arrhythmia initiation requires confirmation in prospective studies. 

By design, our modeling and simulation was relatively sim- 
ple and was limited to a one-dimensional system. This did not 
allow examination of wavebreak, arrhythmia development or 
re-entry, and observations are limited to assumptions of homoge- 
nous properties. The APD and CV restitution models used in 
this paper are in the form of exponential equations, this assumed 
behavior is based on the tendency of experimental data and on 
equations used in previous studies. Biological behavior would 
be expected to deviate from these simple exponential models, 
which in turn will result in discrepancies between simulation 
and clinical measurement. The use of derivations allowing direct 
fitting of ionic models to clinical results, perhaps by using an 



iterative approach, will improve the fidelity of such patient- 
derived simulations. 

The success of our model at predicting local ATs following an 
S3 stimulus makes it tempting to assume that it will retain pre- 
dictive accuracy for S4, S5 and so on. This is unlikely. Though 
the principles demonstrated remain valid, our model did not 
include myocardial memory effects, post-repolarization refrac- 
toriness or anisotropy. These may have increasing importance 
with longer extrastimuli sequences, and ARI restitution dynam- 
ics may change significantly as a consequence of these (Kobayashi 
et al, 1992; Shimizu et al., 2000). Despite a very high agreement 
in predicting S3 ARI restitution obtained from a complex ionic 
model, even these ionic models may still not fully reproduce these 
less studied phenomena. However, our model's simplicity made it 
accessible to the incorporation of clinical data, and the robust pre- 
diction and demonstration of S3 wavefront dynamics (which do 
not require knowledge of S3 ARI restitution) remains impressive. 

FUTURE STUDIES 

This study provides important information to assist electrophys- 
iological characterization of myocardial diseases and phenotypes. 
Some of the techniques we describe may be applicable to clinical 
practice e.g., determining an individual patient heart's boundary 
conditions for functional block as a method for arrhythmic risk 
stratification. All data was acquired in patients with normal hearts 
under baseline conditions. Further studies should examine the 
modulation of changes in conduction/repolarization interactions 
under conditions of autonomic stress or drug interactions. More 
complicated models, such as those based on ionic cellular mod- 
els could be applied to the data, possibly using our method as an 
intermediary step to allow parameter fitting. This could produce 
more detailed and realistic simulations, extending the described 
parameters into 3 dimensions. Furthermore, new technologies 
such as non-invasive body surface mapping may allow acquisition 
of whole heart activation and repolarization dynamic data. 

Specific future uses of the described method may have an early 
application in patients at risk of inherited ventricular arrhyth- 
mias, particularly Brugada syndrome or arrhythmo genie right 
ventricular cardiomyopathy (ARVC). Both of these conditions 
have well-described changes in invasive electrophysiology, yet risk 
stratification is incomplete. Patient specific models could give 
valuable diagnostic information and help quantify the risk of 
arrhythmia in such patients. Pharmaceutical testing may also 
benefit from such analysis, as a method to determine arrhyth- 
mogenicity of subtle changes in electrophysiological responses 
after drug challenge. Furthermore, a bullish use of this method 
would be as an adjunct in early clinical assessments of new 
anti-arrhythmia agents, to lend support to an investigational 
compound's antiarrhythmic properties, or otherwise. 

CONCLUSIONS 

Conduction and repolarization interactions in humans can be 
described in terms of restitution of CV and APD. Simulations 
based on clinically acquired data can be used to successfully pre- 
dict complex activation patterns following sequential extrastim- 
uli. Such modeling techniques may be useful as a method of 
incorporation of clinical data into predictive models. 
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